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The weak version of universality in turbulence refers to the independence of the scaling exponents 
of the nth order strcuture functions from the statistics of the forcing. The strong version includes 
universality of the coefficients of the structure functions in the isotropic sector, once normalized by 
the mean energy flux. We demonstrate that shell models of turbulence exhibit strong universality 
for both forced and decaying turbulence. The exponents and the normalized coefficients are time 
independent in decaying turbulence, forcing independent in forced turbulence, and equal for decaying 
and forced turbulence. We conjecture that this is also the case for Navier-Stokes turbulence. 



I. INTRODUCTION 

The statistical theory of fluid turbulence is concerned 
with correlation functions of the turbulent velocity vec- 
tor field u(r,t) where r is the spatial position and t the 
time Since the velocity field is a vector, multi-point 
and multi-time correlation functions are in general ten- 
sor functions of the vector positions and the scalar times. 
Naturally such functions have rather complicated forms 
which are difficult to measure and to compute. Conse- 
quently, almost from its very beginning, the statistical 
theory of turbulence had been discussed in the context 
of an isotropic and homogeneous model. The notion of 
isotropic turbulence was first introduced by G. I. Taylor 
in 1935 2J. It refers to a turbulent flow, in which the sta- 
tistical averages of every function of the velocity field and 
its derivatives with respect to a particular frame of axes is 
invariant to any rotation in the axes. This is a very effec- 
tive mathematical simplification which, if properly used, 
can drastically reduce the mathematical complexity of 
the theory. For this reason, it was very soon adopted by 
others, such as T. D. Karman and L. Howarth who 
derived the Karman-Howarth equation, and A. N. Kol- 
mogorov |4|,|5j who derived the 4/5 law. In fact, most of 
the theoretical work in turbulence in the past sixty years 
had been limited to the isotropic model. 

Within the homogeneous and isotropic model there de- 
veloped the notion of universality of turbulence. By uni- 
versality, we mean the tendency of different turbulent 
systems to show, for very large Reynolds numbers IZe, 
the same small-scales statistical behavior when the mea- 
surements are done far away from the boundaries. The 
statistical objects (see bellow for definitions) exhibit ap- 
proximately the same scaling exponents whether they are 
measured in the atmospheric boundary layer, in a wind 
tunnel or in a computer simulation, provided they are 
measured far from the boundaries. Moreover, the accu- 
mulated experimental knowledge over the years indicated 
that not only in forced, stationary turbulence, but also 
in decaying turbulence, there is a regime of time where 
the statistical objects exhibit the same scaling properties. 
This phenomenon was explained by the widely sepa- 
rated time scales ("eddy turn over times") that charac- 



terize large and small length scales in turbulence. While 
turbulence was decaying on the time scale of the large ed- 
dies, the small one had ample time to reach an "energy- 
flux equilibrium" that in terms of scaling behavior was 
indistinguishable from forced turbulence. Thus there ex- 
ists a wide-spread belief that at least from the point of 
view of scaling exponents, forced and decaying turbu- 
lence are in the same universality class, sharing the same 
scaling exponents of the corresponding correlation func- 
tions. 

To actually prove this type of universality in experi- 
ments and simulations is however far from straightfor- 
ward. To achieve reasonable precision in the measure- 
ment of scaling exponents one needs large ranges of scales 
where scaling prevails, and this entails large Reynolds 
numbers. Unfortunately large Reynolds numbers are 
available usually when anisotropic effects are large, like 
in the atmospheric boundary layer or in large wind tun- 
nels. Direct Numerical Simulations (DNS) can be used 
to eliminate anisotropy almost completely (up to lattice 
anisotropy which are unavoidable in simulations), but 
they are limited to relatively low IZe, notwithstanding the 
very recent simulations of size 4096 3 ©]. Decaying tur- 
bulence is even harder to characterize precisely, since the 
effective Reynolds number decreases in time. Thus ac- 
tual measurements of scaling properties are fraught with 
difficulties, corrections to scaling, effects of anisotropy 
and what not. As a result, over the years and also 
very recently, it was proposed that structure functions 
in forced and decaying turbulence have different expo- 
nents. In this paper we take a strong stand, proposing 
that the universality that actually exists in turbulence is 
even stronger than what has been anticipated so far. 

To make our point clear, recall that the statistical de- 
scription of fully developed turbulence employs correla- 
tion functions and structure functions. These are ensem- 
ble average of velocity differences across a length-scale R. 
In the theoretical studies of turbulence the two most com- 
mon ensemble averages are over realizations of the forcing 
when one studies forced turbulence, or over initial condi- 
tions when one studies decaying turbulence. The longi- 
tudinal structure functions are the simplest such objects, 
being moments of the longitudinal components of the ve- 
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locity difference between two points. We will denote the 
longitudinal structure functions in forced and decaying 
turbulence by S p and F p respectively, with the precise 
definitions 

Sp (R) = {{[u(r + R,t)-u(r,t)]-*y)t , (1) 
F P (R, t) ee ({[u(r + R,t)- u(r, t)} ■ |}*» ), . (2) 

Here u(r,t) is the velocity field measured at point r at 
time t. (. . .)f and (. . .); stand for ensemble averaging 
over the forcing and the initial conditions respectively. 
In writing Eq. Q we assumed that the forcing is station- 
ary in time, homogeneous and isotropic, and thus S n is 
a function of the scalar R only. In writing Eq. (J2J) wc 
assumed that the initial condition are homogeneous and 
isotropic. Of course the decaying structure functions are 
by definition time dependent. The widely spread belief 
IH H, HI] is for R values in the inertial range of tur- 
bulence (much smaller than the forcing scale but much 
larger than the dissipation scale), the scaling exponents 
Q p that characterize S P (R), i.e. S P (XR) — X'' P S P (R), 
are the same as the scaling exponents that character- 
ize F p (R,t) for a given value of t. Of course also here 
R should be well in the (time dependent) inertial range 
and t should be neither too small nor too large In 
the sequel we refer to the identity of the only scaling 
exponents of these two sets of objects (if it exists) as 
"the weak version of universality" . As mentioned above 
the existence of the weak version of universality is by no 
means accepted by everybody in the field of turbulence. 
Since there is no proof of this universality, doubts of its 
existence linger, and for example in Q it was concluded 
that the scaling exponents of the two families of statis- 
tical objects are not the same. We note however that in 
the same paper it was stated that the scaling exponents 
of the longitudinal and transverse structure functions are 
also not the same. It was shown recently however that 
such statements stem from incomplete treatments of the 
effects of anisotropy [2 El j leaving hope that the weak 
version of universality is still correct. 

In fact, in this paper we will propose that not only the 
weak version of universality is correct, but in fact also 
a "strong version of universality" is applicable. By the 
latter we mean that once properly normalized, the struc- 
ture functions F p and S p agree not only in exponents 
but also in amplitudes. In the context of the 2nd order 
structure function this is not a new statement. The uni- 
versality of (G and C2) was already stated in the 80 's 
by Yaglom [14| and Kader [l5j. Analyzing hundreds of 
experiments made in different flows under different con- 
ditions, Sreenivasan in 1995 came to the conclusion that 
"the Kolmogorov constant C2 is more or less universal, 
essentially independent of the flow as well as the Reynolds 
number (for R\ > 50 or so),., with the average value of 
C2 ~ 0.53 with a standard deviation of about 0.055" [lfij . 
Nevertheless the universality of C2 and C p for p < 4 is 
still under debate. For example, very recently in P"t| 



the authors argued on the basis of a 256 3 DNS "in favor 
of an "exponents only" universality scenario for forced 
turbulence" . 

We believe that this strong version of universality was 
never stated before, and the common thinking is that 
amplitudes depend in a non-universal way on details of 
the forcing or the preparation of the decaying turbulence. 
While true, we will argue that the freedom afforded by 
such details is very limited, amounting at the end to just 
one free number, which, once taken into account, the 
strong version of universality applies. 

Besides being an issue of fundamental importance to 
turbulence, there is another reason for returning at this 
time to the correspondence between force and decaying 
turbulence. The reason is that the riddle of anomalous 
scaling of correlation and structure functions in forced 
turbulent advection (passive and active) had been solved 
recently. First in the context of the non-generic Kraich- 
nan model of passive scalar advection |18| . and then, in 
steps, for passive vectors 0,|2(j , generic passive scalars 
and vectors 0, |2 2| , p j| and finally for generic active 
scalar and vectors Hill! 111. The common thread of 
this advance is that anomalous scaling is discussed in the 
context of the decaying (unforced problem), in which one 
shows that there exist Statistically Preserved Structures 
(eigenf unctions of eigenvalue 1 of the appropriate prop- 
agator of the decaying correlation functions). The de- 
caying problem is independent of forcing, and one shows 
that the statistics of the forced problem is dominated 
by the same Statistically Preserved Structures that are 
identified in the decaying problem. The calculation of 
the anomalous exponents boils down then to calculating 
eigenf unctions of linear operators. In these problems the 
correspondence between the decaying and forced statis- 
tics is proven mathematically or demonstrated beyond 
reasonable doubt by careful numerics. A crucial ingre- 
dient in all this progress is that turbulent advection is 
described by linear pde's. There is therefore an urgent 
question how to translate (if it is possible) the newly ac- 
quired insights to the non-linear turbulent problem itself, 
be it the Navier-Stokes equations or any of the shell mod- 
els that were frequently discussed recently in the context 
of anomalous scaling. In this paper we make a step in this 
direction, analyzing the decay of the Sabra shell model 
[2^ and showing numerically that the statistics of the de- 
caying state and the forced turbulent state are the same 
in exponents and in amplitudes up to one freedom (the 
time dependent mean energy). We opt to work with the 
shell model rather than the Navier-Stokes equations sim- 
ply because the accuracy required for our aims exceeds 
the available scaling ranges and decay times for the lat- 
ter. We express a strong belief that very similar results 
can be demonstrated also for Navier-Stokes turbulence. 
Indeed, in a future publication we will present the the- 
ory that stands behind the present numerical findings and 
demonstrate that the basic structure of that theory is the 
same for shell models and the Navier-Stokes equations. 

The paper is organized as follows: in Sec. [H] we intro- 
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duce the shell model and the numerical simulations that 
we perform. We present the data for the energy decay 
and explain what is the time domain for which we should 
compare the decaying and the forced statistics. In Sec. lIIII 
we present the results for forced structure functions for 
different types of forcing. In determining the exponents 
and the amplitudes of these functions one has to be extra 
careful - we explain that one needs to find fits to functions 
throughout their range of existence. It is not enough to 
plot log-log plots for the inertial range. In Sect. II VI we 
present the data for the decaying correlation functions, 
and explain how to find their exponents and amplitudes 
once the time dependence is taken into account. We ex- 
plain theoretically that the decaying structure functions 
contain sub-leading contributions that decay fast toward 
small scales and do not affect the leading scaling expo- 
nents. The results of our calculations are summarized in 
Table 1 which is the central result of this paper, giving 
strong support to the conjecture of strong universality. 
In Sec. we present a summary and some concluding 
remarks. 



II. MODEL AND ENERGY DECAY 
A. Model and objectives 

The Sabra shell model [27|], like all shell models, is a 
reduced dynamical description of turbulence in terms of 
complex variables u n which represent velocity amplitudes 
associated with wavenumber k n = fcoA™. The equations 
of motion are 

= i(ak n+ iu n+2 u* l+1 + &fc„u„ +:L u*_i (3) 
-c/c„_iu„_iu„_ 2 ) - vk\u n + /„ , 

where the symbol * stands for complex conjugation and 
v is the "viscosity" . The coefficients a, b and c are chosen 
such that a + b + c = 0. This guarantees the conservation 
of the "energy" 



E 



(4) 



in the inviscid (y = 0) forceless limit. As it is well known, 
the Sabra model has a second quadratic invariant, anal- 
ogous to the helicity in fluid mechanics, of the form 



h = 5>/ c y>„| 5 



(5) 



In this paper we will compare the statistics of the 
forced solution (with the forcing /„ restricted to the first 
and second shells, n = 1, 2) to the statistics of the decay- 
ing problem with /„ = for all n. The comparison will be 
presented in terms of the time-independent forced struc- 
ture functions S n and time dependent decaying structure 



functions F n defined as follows: 

S 2 (k n ) = (K| 2 ) f ; F 2 (k n ,t) = (Kl 2 ); , 
S 3 (k„) = Im(w„_iM„u* +1 ) f ; (6) 
F 3 (k n ,t) = Im(w„_iM„u* +1 )i , 

S^K) = (Ki| 4 )f 5 F 4 (k n ,t) = (|w„| 4 )i , 
S 6 (k n ) ee (K| 6 ) f ; F 6 (k n ,t) = (Kl 6 ); . 

Here (. . .)f and (. . .); represent averaging with respects 
to realizations of the forcing and the initial condition 
respectively for the forced and decaying problem. 

The main result of the present work is that in this 
model there exists strong universality. This means that 
in the bulk of the inertial interval the "decaying" struc- 
ture functions F p (k n ,t) take on the form 



ii(t) 



p/3 



\-nCp 



(7) 



with the same anomalous scaling exponents £ p and the 
same dimensionless constants C p , as in the scaling laws 
of the "forced" structure functions 



£f 

h 



IP/3 



X -nC P 



Here 



Ei(t) = (e n (t))i , £f = (s n (t)){ 



(8) 



(9) 



In these equations the instantaneous value of the energy 
flux going through nth shell (in a given realization) is 

e„(i) = 2k n - a \lvt\{u n {t)u n+1 {t)u* n+2 \ (10) 
+ clvc\{u n -i{t)u n {t)u* n+l } . 

(for more details, see [13 )• The only difference between 
Eqs. Q) and ||HJ is that the energy flux e"i(i), averaged 
over the statistics of the initial conditions, is decaying in 
time, while the energy flux E{, averaged over the statis- 
tics of the forcing, is time independent. Eqs. Q) and © 
imply that the probability distribution function (PDF) 
of the velocity fluctuations u n on scales k n within the in- 
ertial range in decaying turbulence can be obtained from 
the corresponding PDF in stationary forced turbulence 
(and vise versa). This is achieved simply by the inter- 
change Si(t) <^> £f. Moreover, strong universality means 
that the PDF dependence on e~i(i) for the fine scales in 
decaying turbulence is independent of the initial condi- 
tions, providing that the Reynolds number IZe 3> 1 and 
all the initial energy is concentrated in the region of large 
scales (first shells). Similarly, the PDF dependence on Ef 
for the fine scales in forced turbulence is independent of 
the statistics of forcing for TZe ^> 1, provided that the 
forcing is concentrated in the region of large scales. 

The rest of this paper is devoted to substantiating 
these (strong) propositions, which can be summarized 
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in the terms of the existence of a probability distribution 
function 

V{ui,...,u N ,t) = [v(t)]~ P{xi, . . . ,x N ) , (11) 

in which = Ui/v(t) and v(t) oc [e„(t)] 1 / 3 is the corre- 
sponding velocity scale. 

B. Simulations 

The calculations presented bellow were carried out for 
the Sabra model with 28 shells, A = 2, a — 1, b = c = 
—0.5, ko = 2~ 4 and v = 10~ 7 . In our simulations we 
employed two different types of forcing. We denote them 
as Forced 1 and Forced 2. Forced 1 has white noise added 
to the equation of the first shell. Forced 2 is forced by 
a Gaussian force on the first shell which is correlated 
exponentially in time. In both cases the amplitude of 
the force in the first shell was chosen /i = 0.01, while 
the forcing amplitude in the second shell was adjusted 
to reduce the helicity input /2 = \J (— c/a)/i, for more 
details, see [27|, pp. 1813 and 1815. 

In the decaying case, the total initial energy E in the 
two first shells was kept constant. The amplitude of the 
first two shell velocities were defined as |mi(0)| 2 = ccEq 
and |u2(0)| 2 = (1 — o)Eq with a random, uniformly dis- 
tributed in the interval [0, 1]. The phases in both shells 
were random, uniformly distributed in the interval [0, 2n] . 
A 4 th order Runge-Kutta scheme with adaptive time step 
was applied. The total energy decay was followed for 9 
decades in time. The statistical objects were accumu- 
lated during 5 decades in time. 

Two recording schemes for the decaying turbulence 
were applied. In one case (denoted as Decay 1), the data 
were recorded starting after a short transient time with 
Eq = 10. For this case, the data was averaged over 13200 
initial conditions. In another case ( Decay 2), the data 
were recorded when the energy in each realization have 
reached the value E = 0.1 with Eq = 5. This data was 
averaged over about 37000 initial conditions. The decay 
of the total energy for two cases, plotted with an appro- 
priate time shift, are shown in Fig. ^ It is clear that the 
two schemes are equivalent for the study of the advanced 
stages of the decay. 

C. The law of energy decay 

We first discuss the total energy decay, where the total 
energy E(t) is defined as 

£W = Eki 2 - ( 12 ) 

n 

In the Navier-Stokes case the law of energy decay 
had been intensively studied following the influential 
works of Taylor El, Kolmogorov 0, Batchelor and 
Townsend [l(J, l28| . For recent development, see, e.g. 
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FIG. 1: The total energy decay from Decay 1 data (dashed 
line) and from Decay2 data ( solid line). The lines coincide 
within the line width. 

Ref. Q and references therein. It was found that E(t) oc 
t~ n with the decay exponent, n, ranging between 1 and 
5/2. This large degree of uncertainty stems from difficul- 
ties in pinpointing the energy spectrum at scales larger 
than the energy containing scale L. It is also not easy 
to determine how L depends on time. Take for exam- 
ple the case of grid turbulence in a wind tunnel. Im- 
mediately behind the grid L is of the order of the mesh 
size. It increases however with the distance from the 
grid. Downstream L may saturate at the wind tunnel di- 
ameter. In this regime the phenomenological analysis Q 
predicts n — 2, which is a number that is not in con- 
tradiction with experiments 11]. The same prediction 
(n = 2) was reached in DNS of the Navier Stokes equa- 
tion and for the GOY shell model [gjj. This prediction 
was shown to be in agrement with numerical simulation 
in which L is time independent due to the special choice 
of initial conditions. 

For the sake of completeness we review the theoret- 
ical analysis of [29|. and show that our simulations of 
the Sabra shell model are in excellent agreement with 
its predictions. Consider a decaying solution with the 
energy initially concentrated, say, in the first two shells, 
see Fig. [21 upper panel. Time is measured in natural 
time units T which are determined by the characteristic 
time of the first shell, T = 1/SfcoV^ (Eo is the total 
initial energy). One sees that during one T the energy 
cascaded down to the 6th shell, and during 2T down to 
the 12th shell. At later times the cascade process acceler- 
ates, and the energy goes from the 13th shell to "infinite" 
shells during a time that is roughly between 2.75 T and 
3 T. As expected, the completion of the cascade process 
requires a finite time T* of a few units T. In our case 
T* « 3T. 

For t > T* the total energy of the system E(t) begins 
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FIG. 2: The behavior of F<z, Development of the energy cas- 
cade ( upper panel). The decay phase ( lower panel). 



to decay with 



dE(t) 
dt 



(13) 



Accordingly, the time dependent Reynolds number 
lZe(t), which is proportional to \J E(t), decreases, and 
the viscous cutoff kd(t) is moving toward smaller shell 
numbers, as is shown in Fig. El lower panel. For example, 
fc d (10 5 T) « 2 12 , fc d (10 8 T) w 2 3 and the inertial interval 
almost disappeares. For larger times all the energy is 
contained in the first shell, and it decays exponentially, 



E(t) oc exp[-2vk\t] 



(14) 



following the linear part of the equation of motion for the 
first shell. 

For intermediate times, which in our simulations span 
the eight orders of magnitude for 3 T < t < 10 8 T, the 
slope of plots of log2i 7 2(fcn) vs n remains more or less con- 
stant. This is a manifestations of the time independence 



FIG. 3: The decay of the total energy E = ^ n \u„\ 2 averaged 
over 3200 initial conditions. The dashed line corresponds to 
the decay law, Eq. with E* = 12.15 and U = 8.35 T. 



of the scaling exponents. Taking this as a fact, one imme- 
diately see from Eq. for p = 2 that E(t) oc [£i(i)] 2/3 
and hence £\{t) oc [E(t)] 3 ^ 2 . Note that this result is inde- 
pendent of the precise value of the scaling exponent £2 , 
anomalous or not. Thus Eq. (|13fl can be presented as 



dE(t) 
dt 



(15) 



with a pre-factor k, which may be expressed via the pa- 
rameters of the shell model, k , (2 an d C%. Approxi- 
mately, k — k{) . In our calculations, both Decay 1 and 
Decay 2, k = 0.0687, while fc = 0.0625. The solution of 
Eq. JTSJ is 



E(t) = E, 



it 



(16) 



where = 2/ ' k-\/E* and E„ is the integration constant. 
The results of the numerical simulations for the total en- 
ergy, E(t), (cf. Fig. Ill CI solid line), are in excellent 
agreement with Eq. (|16fl . which is shown in the figure 
as a dashed line. The total energy decay was followed 
for 9 decades in time. After 1.5 decades of transient 
behavior, the decay of total energy follows very closely 
the t~ 2 law, until at about 6 decades the viscous scale 
reaches the first shells and the decay become exponential 
in agreement with Eq. Q14[l. 



III. FORCED STRUCTURE FUNCTIONS 

In this section we present results for the forced struc- 
ture functions. As far as the scaling exponents are con- 
cerned, there is not much novelty in this section, the 
exponents are basically the same as those reported in 
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FIG. 4: The structure function for two types of the forc- 
ing. Large symbols, shown in the figure correspond to Forced 
2 data. The dashed lines are the fits for the corresponding 
structure functions. Small black dots denote Forced 1 data. 
The fits for these data are not shown. Both sets were normal- 
ized by their respective total energy. 



the shells taken into account except the first two and the 
last three shells, to minimize the boundary effects. 

The quality of the fit may be seen in Fig.0J The forced 
data are shown normalized by the respective total energy, 
but not compensated by k^ . Then the different structure 
functions are separated and data for Force 1 and Force 2 
cases may be distinguished. 

The fit procedure allows us to express the structure 
function in the inertial range as S p = A p EP^ki P . To 
calculate the coefficients C p we have now estimate the 
value of the energy flux £f[ see Eq. (JS}]. We use the 
exact result for £3 to express E{ via A 3 and parameters 
of the model Eq. (|2(J[1 . The coefficient C p of the structure 
functions, other than S3, may be therefore written as: 

Cp = [2a 3 (a- C )K3 ■ (19) 

The results are summarized in the Table [I] Before we 
discuss this Table, which is the central result of this pa- 
per, we turn to the analysis of the decaying structure 
functions and add their analogous results to the Table as 
well. 



a number of previous publications. The new aspect 
stressed here is that of the coefficients C p of the struc- 
ture functions, cf. Eq. (JSJ. We demonstrate that these 
coefficients are independent of the forcing, and that the 
scaling form proposed in Eq. © is indeed universal. 

To get an accurate determination of the scaling expo- 
nents and of the coefficients, and to be able to demon- 
strate strong universality, it is mandatory to fit the mea- 
sured data to model functions that contain the presumed 
dissipative behavior. Failing to do so results in inaccu- 
racies that may lead to infinite confusions. As a first 
step in analysis of the data we fit the normalized struc- 
ture functions S p /E p l 2 with the fit formula Eq. (44) from 



V k d , p J 



(17) 



where 



^p? Cp 7 ^p j Mp s ^d,p 5 



or 



nd,p = log 2 k d , P 



are the fit parameters and x = log A (l + \/5)/2 is the 
exponent of the viscous range. The parameters A p ,£ p 
determine the behavior of P p (k n ) in the inertial inter- 
val, fcd,p determines the viscous cutoff. The "auxiliary 
parameters" a p ,/z p correct the behavior in the transient 
inertial-viscous region. To obtain the best fit we mini- 
mize the error function: 



£ 



ioglO Pp( k r> 

logic S p (k n 



(18) 



where S p refers to the numerically obtained data. Both 
sets of forced structure functions data were fit with all 



IV. DECAYING STRUCTURE FUNCTIONS 

In this section we present results for the decaying struc- 
ture functions, including the numerical support for the 
strong universality proclaimed in Eq. J7J. We caution 
the reader (and whoever wants to repeat these calcula- 
tions in other systems, including Navier-Stokes decaying 
turbulence), that the issue is fraught with sub-leading 
contributions, even in the isotropic sector. One obvious 
sub-leading term is provided by the rate of change of 
Fp(kn,t) which is coupled, via the infinite hierarchy of 
equations, to terms involving F p+ i(k m ,t), with m of the 
order of n. To sec this phenomenon clearly and to learn 
how to take it into account we discuss first the case of 
the 3rd order structure function which can be dealt with 
analytically. 



A. Sub-leading corrections to the scaling of 
decaying turbulence 

The easiest case for theoretical analysis of the scaling 
behavior of the decaying structure functions F p (k n ,t) is 
the case p = 3. for which in the forced case Ss(k n ) is 
known exactly [2Jj. For simplicity we will discuss here 
the helicity-free case, for which 



S 3 (k n ) = £ f /fc„(c- a) 



(20) 



We will show now that in the decaying case strong univer- 
sality is realized, but only well within the inertial range. 
In the vicinity of the energy containing scales there are 
significant sub-leading corrections caused by the time de- 
pendence of £i(t). 
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TABLE I: Universal coefficients C p and scaling exponents measured from "the best fit" on numerical data. The auxiliary fit 
parameters in Eq. 1171 are found to be in the intervals a p ~ 0.5 -f- 2 and /j p ~ 0.6 -f- 2.3. The error bars for each parameter 
correspond to the error function £ , Eq. (1181 ) equal ~ v2^min with all other parameters set to their optimal values. 
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1.72 ±0.01 


15.87 


Forced 2 




0.73 ±0.07 


0.728 ± 0.006 


17.0 


0.60 ±0.08 


1.254 ±0.008 


16.24 


0.61 ±0.15 


1.72 ±0.01 


15.87 


Decay 1 


20 


0.72 ± 0.05 


0.728 ±0.006 


17.2 


0.62 ±0.06 


1.254 ± 0.008 


17.0 


0.79 ±0.15 


1.72 ±0.01 


16.4 




100 


0.72 ± 0.06 


0.728 ± 0.006 


16.2 


0.62 ±0.07 


1.254 ±0.008 


16.0 


0.79 ±0.15 


1.72 ±0.01 


15.0 




10 3 


0.72 ±0.07 


0.728 ± 0.006 


13.5 


0.61 ±0.08 


1.255 ±0.008 


13.2 


0.77 ±0.15 


1.72 ±0.01 


13.0 




10 4 


0.72 ± 0.08 


0.728 ± 0.006 


11.0 


0.61 ±0.09 


1.256 ±0.008 


10.2 


0.75 ±0.15 


1.73 ±0.01 


9.8 


Decay 2 


100 


0.73 ±0.1 


0.73 ±0.01 


15.5 


0.6 ±0.2 


1.25 ±0.01 


15.2 


0.70 ±0.20 


1.72 ±0.02 


14.3 




10 3 


0.74 ±0.1 


0.73 ±0.01 


13.3 


0.6 ±0.15 


1.26 ±0.01 


13.1 


0.66 ±0.20 


1.72 ±0.02 


12.4 




10 4 


0.72 ±0.1 


0.73 ±0.01 


10.7 


0.6 ±0.2 


1.25 ±0.01 


10.5 


0.72 ±0.25 


1.73 ±0.03 


9.8 



To find these corrections, consider the equation of mo- 
tion of the 2'nd order structure in the inertial interval 
(i.e. for k n < k d ) {Eq. (9) of Ref. 

dF ^ kn : t] = a\F s (k n+1 ,t) + bF 3 (k n ,t) + ^F 3 (fc n _ l7 <) . 
lk n dt X 

(21) 

In the stationary case, the LHS of this equation vanishes 
and Eq. (|2UI) is a solution. In the decaying case, let 

F 3 (k n ) = F 3 {a) (k n ,t) + SF 3 (k n ,t), (22) 

*f(A») = • 

k n {c - a) 

In the intermediate time regime F2(k n , t) cx (t + i*)~ 2 , 

dF 2 (k„,t)/2k„dt^ -F 2 {k n ,t)/k n {t + t,) . (23) 

Comparing with Eq. I|21|l . we see that the leading solu- 
tion for F 3 gains a sub-leading term SF 3 (k n , t) cx kn^^ 
or more precisely 



5F 3 (k n ,t) 



F 2 (k n ,t) 



(t + U) k n (a + b + c\&) 



(24) 



Notice, that F^ (k n ,t) and SF 3 (k n , t) have the same time 
dependence, cx (t ± i*)~ 3 , but different scaling. As a 
result, their ratio is time independent; the sub-leading 
term does not become relatively smaller in time. On the 
other hand it decays relatively to the leading term as k n 
increases: 



SF 3 (k n ,t) 1 
—) ° C_ 



F^(k n ,t) 



\nC 2 



(25) 



Although the 3rd order structure function is easiest to 
handle, it is clear that there will always be a subleading 
term added to F p from the time derivative of F p -% which 



appears in the infinite hierarchy of equations. Since these 
equations always have k n on the RHS, one can immedi- 
ately guess the general form of the correction to scaling 
for F p , i.e. 



1 



SF p (k n ,t) 

£__ _ ________ 

F p (k n ,t) \n(l+<P-i-< P ) 



(26) 



This means that the sub-leading term of the p-order 
structure functions decreases toward small scales roughly 
as A~ 2n / 3 for "normal" K41 scaling, and somewhat slower 
for anomalous scaling. Strong universality of the turbu- 
lent statistics is thus expected only deeply in the inertial 
interval. 



B. Numerical results 

All calculated statistical objects were normalized by 
the total energy, Fp{k ni t)/EP/ 2 (t). All the normal- 
ized, compensated decaying structure functions show a 
plateau. To enrich the statistics the data were first nor- 
malized by the total energy and then averaged over one 
tenth of a temporal decade. For the same reasons that 
were explained in the case of the forced objects, the fit 
region for the decaying structure functions was chosen 
from n — 3 to n = n^.p + 5 for the Decay 1 data and from 
ri = 5 to n = rid,p ± 5 for the Decay 2. For t — 10 5 T 
rid ~ 7 and only very few shells may be considered as the 
"inertial interval" . Therefore the fit parameters become 
unreliable. The quality of the fit can be seen in Fig. _for 
t between 20 T and 10 4 T. For t<20T the flux equilib- 
rium cannot be guaranteed and the coefficients C p and 
the exponents may be not universal. This is definitely 
the case for t < 10 T. The decaying structure functions 
are plotted normalized by the total energy and compen- 
sated to emphasis the fact that main effect is the shift of 



J 



8 




-5 

UJ 

w M -15- 



«- -9-* -«-« -*-«-«-«-«- =9- 9 -fc g ^ 

^ \ ^ 



CD 



-20- 
-25 
-30- 
-35 



o 20 

. 100 

« 10 3 

< 10 4 



10 



V? 



15 



20 




-5 

r -10 
ii 

>>-15 

x 

£ c -20 

CM 

3 -25 
-30 
-35 



o 20 

* 100 

<> 10 3 

< 10 4 



' •« -4- -•- ♦ H»- ♦ -f • -•.£>.. 



I 1 



1 q 



10 15 
n 



20 



-5 

-10- 



« -15- 



C\J 



-20- 

-25 

-30- 

-35^ 




o 20 

* 100 

o 10 3 

< 10 4 



10 



15 



20 



x 

^20- 

O) 

o 
_1 

-25 



-35 b 




-5 
10 

15- 



o 


20 




100 


<> 


10 3 


> 


10 4 



10 15 
n 



x * \o 



1 I 1 



f 0+ 



20 



FIG. 5: The normalized compensated decaying structure functions for Decay 1 data averaged over 13200 initial conditions and 
over one tenth of the decade in time. The symbols show the calculated data at different times( defined in the legend). The 
dashed lines denote the fits for the corresponding structure functions. The best fit scaling exponents, used for the compensation 
(£2 = 0.728, = 1-, C4 = 1.254 and £e = 1.72) are the same as in the forced case. The Decay 2 data show similar behavior. 



V. SUMMARY AND DISCUSSION 

In summary, we analyzed, using the Sabra shell model 
of turbulence, "forced" structure functions S p (k n ) for 
two types of forcing, and "decaying" structure func- 
tions F p {k n , t) for two types of initial conditions, (called 
Decay 1 and Decay 2). For Decay 1 we considered 
four times, which differ in order of magnitudes (t = 
20 T, 10 2 T, 10 3 T, and 10 4 T, where T is the character- 
istic time of the 1st shell). For Decay 2 we considered 
three different times t = 10 2 T, 10 3 T, and 10 4 T. In all 
these cases we found the scaling exponents, £ p , and the 
dimensionless amplitudes, C p , for the three even orders 
p = 2, 4 and 6. The results are collected in Table 1 to- 
gether with our estimates of the error bars. 

We can state that our results support the conjecture of 



strong universality within the numerical accuracy. Con- 
cerning the fact that before the data analysis presented 
above the raw results contained objects differing by or- 
ders of magnitude, the degree of precision of the identity 
of the amplitudes C p , and the exponents Q p , shown in Ta- 
ble 1 should be taken very seriously. We propose that all 
the results presented by previous authors with negative 
indications about universality (even of the weak type) 
stem from problems in handling the corrections to scal- 
ing, either from anisotropy or from dissipative or other 
boundary effects. 

It should be stressed at this point that the strong 
universality observed here is not expected in the much 
simpler problem of turbulent advection. The difference 
stems from the linearity of the advection problem vs. 
the nonlinearity of the Navier-Stokes problem and its 



9 



shell counterpart. In the linear advection problem one 
finds equations for the statistical objects that decouple 
for each order p. The independent p-order equations de- 
termine the anomalous scaling exponent £ p from solvabil- 
ity conditions, leaving the amplitude C p to be found by 
matching the scale invariant correlation function in the 
inertial interval with its non-universal "boundary condi- 
tions" at energy contained scales. The amplitudes C p de- 
pend therefore on the details of the non-universal forcing, 
and the statistics of the turbulent advection problem may 
exhibit weak universality only. In contrast, the nonlinear- 
ity of the Navier Stokes equations and their shell-model 
counterparts leads to coupled, hierarchical equations for 
all the p-order statistical objects that have to be solved 
for simultaneously. This rigid structure allows much less 
freedom than the linear advection problem, leading to 
the possibility of strong universality. 

Finally, we comment on a possible theoretical sup- 
port for the strong universality conjecture. We propose 
that a necessary condition for strong universality is the 
locali ty o f interaction, which allows to formulate (see 
Refs. |3(j,|3l],|32) the hierarchy of equations in terms of 
inertial-range objects only. The locality of energy trans- 



fer over scales, which is built in the shell models of tur- 
bulence, is an assumption in the Richardson-Kolmogorov 
cascade picture of turbulence, see, e.g. 0, 0, E|. T he lo- 
cality of interaction was demonstrated in Ref . 33] , using 
the Belinicher-L'vov transformation of the Navier-Stokes 
equations |34j , which allows to eliminate from the theory 
the sweeping effect. Once we have a theory in terms of 
inertial-range objects, it is quite acceptable that ampli- 
tudes should be universal as well, up to an overall single 
parameter which is the energy flux per unit time and 
mass. An elaboration of these ideas will be presented in 
a future publication. At this point we finish with the 
conjecture that strong universality is a property shared 
also by Navier-Stokes turbulence. 
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